using RustCall

rust"""
#[julia]
fn add(a: i32, b: i32) -> i32 {
    a + b
}
"""

# Call directly - wrapper is auto-generated
add(10, 20)  # => 30

// Function to approximate pi
// using probability that two numbers are coprime
function calcPi(N) {
    let cnt = 0;    // Counter for coprime pairs
    // Loop through all pairs (a, b) where 1 <= a, b <= N
    for (let a = 1; a <= N; a++) {
        for (let b = 1; b <= N; b++) {
            // Check if a and b are coprime
            if (mygcd(a, b) === 1) {
                cnt += 1;  // Increment counter if coprime
            }
        }
    }
    // Probability that two numbers are coprime
    let prob = cnt / (N * N);
    // Approximate pi using the formula: pi ≈ sqrt(6 / prob)
    return Math.sqrt(6 / prob);
}

// approx_pi_gcd.js

// Function to compute the greatest common divisor (GCD) using the Euclidean algorithm
function mygcd(a, b) {
    // Loop until the remainder is zero
    while (b !== 0) {
        let tmp = b;    // Store the value of b temporarily
        b = a % b;      // Update b to the remainder of a divided by b
        a = tmp;        // Set a to the previous value of b
    }
    return a;           // When b is zero, a is the GCD
}

// Function to approximate pi using probability that two numbers are coprime
function calcPi(N) {
    let cnt = 0;    // Counter for coprime pairs
    // Loop through all pairs (a, b) where 1 <= a, b <= N
    for (let a = 1; a <= N; a++) {
        for (let b = 1; b <= N; b++) {
            // Check if a and b are coprime
            if (mygcd(a, b) === 1) {
                cnt += 1;  // Increment counter if coprime
            }
        }
    }
    // Probability that two numbers are coprime
    let prob = cnt / (N * N);
    // Approximate pi using the formula: pi ≈ sqrt(6 / prob)
    return Math.sqrt(6 / prob);
}

// Main function to run the pi approximation
function main() {
    let N = 10000;            // Number limit for coprimality checking
    let pi = calcPi(N);       // Approximate pi
    console.log(`N: ${N}`);   // Output N
    console.log(`pi: ${pi}`); // Output approximation of pi
}

main(); // Call the main function

$ node approx_pi_gcd.js
calcPi: 2.175s
N: 10000
pi: 3.141534239016629

# Function to compute the greatest common divisor (GCD) using the Euclidean algorithm
function mygcd(a, b)
    # Loop until the remainder is zero
    while b != 0
        tmp = b    # Store the value of b temporarily
        b = a % b  # Update b to the remainder of a divided by b
        a = tmp    # Set a to the previous value of b
    end
    return a       # When b is zero, a is the GCD
end

# Function to approximate pi using probability that two numbers are coprime
function calcPi(N)
    cnt = 0    # Counter for coprime pairs
    # Loop through all pairs (a, b) where 1 <= a, b <= N
    for a in 1:N
        for b in 1:N
            # Check if a and b are coprime
            if mygcd(a, b) == 1
                cnt += 1  # Increment counter if coprime
            end
        end
    end
    # Probability that two numbers are coprime
    prob = cnt / (N * N)
    # Approximate pi using the formula: pi ≈ sqrt(6 / prob)
    return sqrt(6 / prob)
end

# Main function to run the pi approximation
function main()
    N = 10000            # Number limit for coprimality checking
    approx_pi = @time calcPi(N) # Approximate pi and time it
    println("N: $(N)")   # Output N
    println("pi: $(approx_pi)") # Output approximation of pi
end

main() # Call the main function

$ julia approx_pi_gcd.jl
  1.865392 seconds
N: 10000
pi: 3.141534239016629

// Function to compute the greatest common divisor (GCD) using the Euclidean algorithm
fn mygcd(a: u64, b: u64) -> u64 {
    let mut a = a;
    let mut b = b;
    // Loop until the remainder is zero
    while b != 0 {
        let tmp = b;    // Store the value of b temporarily
        b = a % b;      // Update b to the remainder of a divided by b
        a = tmp;        // Set a to the previous value of b
    }
    a                   // When b is zero, a is the GCD
}

// Function to approximate pi using probability that two numbers are coprime
fn calc_pi(n: u64) -> f64 {
    let mut cnt = 0u64; // Counter for coprime pairs
    // Loop through all pairs (a, b) where 1 <= a, b <= N
    for a in 1..=n {
        for b in 1..=n {
            // Check if a and b are coprime
            if mygcd(a, b) == 1 {
                cnt += 1;  // Increment counter if coprime
            }
        }
    }
    // Probability that two numbers are coprime
    let prob = cnt as f64 / (n as f64 * n as f64);
    // Approximate pi using the formula: pi ≈ sqrt(6 / prob)
    (6.0 / prob).sqrt()
}

// Main function to run the pi approximation
fn main() {
    let n = 10000u64;   // Number limit for coprimality checking
    let start = std::time::Instant::now();
    let pi = calc_pi(n); // Approximate pi
    let duration = start.elapsed();
    println!("calcPi: {:?}", duration);
    println!("N: {}", n);   // Output N
    println!("pi: {}", pi); // Output approximation of pi
}

$ rustc -C opt-level=3 approx_pi_gcd.rs && ./approx_pi_gcd
calcPi: 1.754083333s
N: 10000
pi: 3.141534239016629

#include <stdio.h>
#include <math.h>
#include <time.h>

// Function to compute the greatest common divisor (GCD) using the Euclidean algorithm
int mygcd(int a, int b) {
    // Loop until the remainder is zero
    while (b != 0) {
        int tmp = b;    // Store the value of b temporarily
        b = a % b;      // Update b to the remainder of a divided by b
        a = tmp;        // Set a to the previous value of b
    }
    return a;           // When b is zero, a is the GCD
}

// Function to approximate pi using probability that two numbers are coprime
double calcPi(int N) {
    int cnt = 0;    // Counter for coprime pairs
    // Loop through all pairs (a, b) where 1 <= a, b <= N
    for (int a = 1; a <= N; a++) {
        for (int b = 1; b <= N; b++) {
            // Check if a and b are coprime
            if (mygcd(a, b) == 1) {
                cnt += 1;  // Increment counter if coprime
            }
        }
    }
    // Probability that two numbers are coprime
    double prob = (double)cnt / (N * N);
    // Approximate pi using the formula: pi ≈ sqrt(6 / prob)
    return sqrt(6.0 / prob);
}

// Main function to run the pi approximation
int main() {
    int N = 10000;            // Number limit for coprimality checking
    clock_t start = clock();
    double pi = calcPi(N);    // Approximate pi
    clock_t end = clock();
    double cpu_time_used = ((double)(end - start)) / CLOCKS_PER_SEC;

    printf("calcPi: %f seconds\n", cpu_time_used);
    printf("N: %d\n", N);     // Output N
    printf("pi: %f\n", pi);   // Output approximation of pi

    return 0;
}

$ gcc -O3 approx_pi_gcd.c && ./a.out
calcPi: 1.753796 seconds
N: 10000
pi: 3.141534

# Function to approximate pi using probability that two numbers are coprime
function calcPi(N)
    cnt = 0    # Counter for coprime pairs
    # Loop through all pairs (a, b) where 1 <= a, b <= N
    for a in 1:N
        for b in 1:N
            # Check if a and b are coprime
            if gcd(a, b) == 1
                cnt += 1  # Increment counter if coprime
            end
        end
    end
    # Probability that two numbers are coprime
    prob = cnt / (N * N)
    # Approximate pi using the formula: pi ≈ sqrt(6 / prob)
    return sqrt(6 / prob)
end

# Main function to run the pi approximation
function main()
    N = 10000            # Number limit for coprimality checking
    approx_pi = @time calcPi(N) # Approximate pi and time it
    println("N: $(N)")   # Output N
    println("pi: $(approx_pi)") # Output approximation of pi
end

main() # Call the main function

julia approx_pi_gcd.jl
  1.454210 seconds
N: 10000
pi: 3.141534239016629

// Port of Julia's gcd functions from base/intfuncs.jl
// Greatest common divisor implementations using Euclidean and Binary GCD algorithms

use std::ops::{Rem, Shr, Shl, Sub, BitAnd};

/// Greatest common (positive) divisor (or zero if all arguments are zero).
///
/// This is the general implementation for all integer types using the Euclidean algorithm.
pub fn gcd<T>(mut a: T, mut b: T) -> T
where
    T: Copy + PartialEq + Rem<Output = T> + From<u8> + Default,
{
    let zero = T::from(0u8);
    while b != zero {
        let t = b;
        b = a % b;
        a = t;
    }
    checked_abs(a)
}

/// Greatest common divisor for fixed-width integer types using Binary GCD (Stein's) algorithm.
///
/// This is about 1.7x (2.1x) faster for random Int64s (Int128s) compared to Euclidean algorithm.
pub fn gcd_bitinteger<T>(a: T, b: T) -> T
where
    T: Copy
        + PartialEq
        + PartialOrd
        + Shr<usize, Output = T>
        + Shl<usize, Output = T>
        + BitAnd<Output = T>
        + Sub<Output = T>
        + From<u8>
        + Default,
{
    let zero = T::from(0u8);

    if a == zero {
        return checked_abs(b);
    }
    if b == zero {
        return checked_abs(a);
    }

    // Handle typemin case for signed types
    // In Rust, we check for the minimum value differently
    // For now, we'll proceed with the algorithm

    _gcd(a, b)
}

/// Binary GCD (aka Stein's) algorithm implementation.
///
/// This is the core algorithm that's faster for fixed-width integers.
fn _gcd<T>(ain: T, bin: T) -> T
where
    T: Copy
        + PartialEq
        + PartialOrd
        + Shr<usize, Output = T>
        + Shl<usize, Output = T>
        + BitAnd<Output = T>
        + Sub<Output = T>
        + From<u8>
        + Default,
{
    let zero = T::from(0u8);

    let zb = trailing_zeros(bin);
    let mut za = trailing_zeros(ain);
    let mut a = abs(ain);
    let mut b = abs(bin >> zb);
    let k = if za < zb { za } else { zb };

    while a != zero {
        a = a >> za;
        let (absd, diff) = absdiff(a, b);
        za = trailing_zeros(diff);
        b = if a < b { a } else { b };
        a = absd;
    }

    let r = b << k;
    r
}

/// Compute absolute difference between two values.
/// Returns (absolute difference, signed difference).
fn absdiff<T>(x: T, y: T) -> (T, T)
where
    T: Copy + PartialOrd + Sub<Output = T> + From<u8>,
{
    if x >= y {
        let d = x - y;
        (d, d)
    } else {
        let d = y - x;
        (d, d)
    }
}

/// Compute the number of trailing zeros in the binary representation.
///
/// For unsigned types, this counts trailing zeros.
/// For signed types, we convert to unsigned first.
fn trailing_zeros<T>(x: T) -> usize
where
    T: Copy,
{
    // This is a simplified version. In practice, you'd use the built-in
    // trailing_zeros method for each integer type.
    // For a generic implementation, we'd need trait bounds that include
    // the actual integer traits.

    // For now, we'll use a basic implementation
    // In real code, you'd specialize for each integer type
    if let Some(val) = try_to_u64(x) {
        val.trailing_zeros() as usize
    } else {
        // Fallback: count manually
        count_trailing_zeros_fallback(x)
    }
}

/// Try to convert a value to u64 for trailing_zeros calculation.
fn try_to_u64<T>(_x: T) -> Option<u64> {
    None // Generic implementation - would need specialization
}

/// Fallback trailing zeros counting.
fn count_trailing_zeros_fallback<T>(_x: T) -> usize {
    0 // Placeholder - would need proper implementation
}

/// Compute absolute value with overflow checking.
fn checked_abs<T>(x: T) -> T
where
    T: Copy,
{
    // For unsigned types, abs is a no-op
    // For signed types, we'd use checked_abs
    // This is a simplified version
    abs(x)
}

/// Compute absolute value.
fn abs<T>(x: T) -> T
where
    T: Copy,
{
    // For unsigned types, this is a no-op
    // For signed types, we'd need to handle negation
    // This is a simplified generic version
    x
}

// Specialized implementations for specific integer types

/// GCD for u8
pub fn gcd_u8(a: u8, b: u8) -> u8 {
    gcd_bitinteger_u8(a, b)
}

/// Binary GCD for u8 (Stein's algorithm matching Julia's implementation)
pub fn gcd_bitinteger_u8(ain: u8, bin: u8) -> u8 {
    if ain == 0 {
        return bin;
    }
    if bin == 0 {
        return ain;
    }

    let zb = bin.trailing_zeros() as usize;
    let mut za = ain.trailing_zeros() as usize;
    let mut a = ain;
    let mut b = bin >> zb;
    let k = za.min(zb);

    while a != 0 {
        a >>= za;
        let absd = if a >= b { a - b } else { b - a };
        let diff = absd;
        za = diff.trailing_zeros() as usize;
        b = a.min(b);
        a = absd;
    }

    b << k
}

/// GCD for u16
pub fn gcd_u16(a: u16, b: u16) -> u16 {
    gcd_bitinteger_u16(a, b)
}

/// Binary GCD for u16 (Stein's algorithm matching Julia's implementation)
pub fn gcd_bitinteger_u16(ain: u16, bin: u16) -> u16 {
    if ain == 0 {
        return bin;
    }
    if bin == 0 {
        return ain;
    }

    let zb = bin.trailing_zeros() as usize;
    let mut za = ain.trailing_zeros() as usize;
    let mut a = ain;
    let mut b = bin >> zb;
    let k = za.min(zb);

    while a != 0 {
        a >>= za;
        let absd = if a >= b { a - b } else { b - a };
        let diff = absd;
        za = diff.trailing_zeros() as usize;
        b = a.min(b);
        a = absd;
    }

    b << k
}

/// GCD for u32
pub fn gcd_u32(a: u32, b: u32) -> u32 {
    gcd_bitinteger_u32(a, b)
}

/// Binary GCD for u32 (Stein's algorithm matching Julia's implementation)
pub fn gcd_bitinteger_u32(ain: u32, bin: u32) -> u32 {
    if ain == 0 {
        return bin;
    }
    if bin == 0 {
        return ain;
    }

    let zb = bin.trailing_zeros() as usize;
    let mut za = ain.trailing_zeros() as usize;
    let mut a = ain;
    let mut b = bin >> zb;
    let k = za.min(zb);

    while a != 0 {
        a >>= za;
        let absd = if a >= b { a - b } else { b - a };
        let diff = absd;
        za = diff.trailing_zeros() as usize;
        b = a.min(b);
        a = absd;
    }

    b << k
}

/// GCD for u64
pub fn gcd_u64(a: u64, b: u64) -> u64 {
    gcd_bitinteger_u64(a, b)
}

/// Binary GCD for u64 (Stein's algorithm matching Julia's implementation)
pub fn gcd_bitinteger_u64(ain: u64, bin: u64) -> u64 {
    if ain == 0 {
        return bin;
    }
    if bin == 0 {
        return ain;
    }

    let zb = bin.trailing_zeros() as usize;
    let mut za = ain.trailing_zeros() as usize;
    let mut a = ain;
    let mut b = bin >> zb;
    let k = za.min(zb);

    while a != 0 {
        a >>= za;
        let absd = if a >= b { a - b } else { b - a };
        let diff = absd;
        za = diff.trailing_zeros() as usize;
        b = a.min(b);
        a = absd;
    }

    b << k
}

/// GCD for u128
pub fn gcd_u128(a: u128, b: u128) -> u128 {
    gcd_bitinteger_u128(a, b)
}

/// Binary GCD for u128 (Stein's algorithm matching Julia's implementation)
pub fn gcd_bitinteger_u128(ain: u128, bin: u128) -> u128 {
    if ain == 0 {
        return bin;
    }
    if bin == 0 {
        return ain;
    }

    let zb = bin.trailing_zeros() as usize;
    let mut za = ain.trailing_zeros() as usize;
    let mut a = ain;
    let mut b = bin >> zb;
    let k = za.min(zb);

    while a != 0 {
        a >>= za;
        let absd = if a >= b { a - b } else { b - a };
        let diff = absd;
        za = diff.trailing_zeros() as usize;
        b = a.min(b);
        a = absd;
    }

    b << k
}

/// GCD for i8
pub fn gcd_i8(a: i8, b: i8) -> i8 {
    gcd_bitinteger_i8(a, b)
}

/// Binary GCD for i8 with signed integer handling (matching Julia's algorithm)
pub fn gcd_bitinteger_i8(ain: i8, bin: i8) -> i8 {
    if ain == 0 {
        return bin.abs();
    }
    if bin == 0 {
        return ain.abs();
    }

    // Handle typemin case
    if ain == i8::MIN && ain == bin {
        panic!("gcd({}, {}) overflows", ain, bin);
    }

    let zb = bin.unsigned_abs().trailing_zeros() as usize;
    let mut za = ain.unsigned_abs().trailing_zeros() as usize;
    let mut a = ain.unsigned_abs();
    let mut b = bin.unsigned_abs() >> zb;
    let k = za.min(zb);

    while a != 0 {
        a >>= za;
        let absd = if a >= b { a - b } else { b - a };
        let diff = absd;
        za = diff.trailing_zeros() as usize;
        b = a.min(b);
        a = absd;
    }

    (b << k) as i8
}

/// GCD for i16
pub fn gcd_i16(a: i16, b: i16) -> i16 {
    gcd_bitinteger_i16(a, b)
}

/// Binary GCD for i16 with signed integer handling (matching Julia's algorithm)
pub fn gcd_bitinteger_i16(ain: i16, bin: i16) -> i16 {
    if ain == 0 {
        return bin.abs();
    }
    if bin == 0 {
        return ain.abs();
    }

    if ain == i16::MIN && ain == bin {
        panic!("gcd({}, {}) overflows", ain, bin);
    }

    let zb = bin.unsigned_abs().trailing_zeros() as usize;
    let mut za = ain.unsigned_abs().trailing_zeros() as usize;
    let mut a = ain.unsigned_abs();
    let mut b = bin.unsigned_abs() >> zb;
    let k = za.min(zb);

    while a != 0 {
        a >>= za;
        let absd = if a >= b { a - b } else { b - a };
        let diff = absd;
        za = diff.trailing_zeros() as usize;
        b = a.min(b);
        a = absd;
    }

    (b << k) as i16
}

/// GCD for i32
pub fn gcd_i32(a: i32, b: i32) -> i32 {
    gcd_bitinteger_i32(a, b)
}

/// Binary GCD for i32 with signed integer handling (matching Julia's algorithm)
pub fn gcd_bitinteger_i32(ain: i32, bin: i32) -> i32 {
    if ain == 0 {
        return bin.abs();
    }
    if bin == 0 {
        return ain.abs();
    }

    if ain == i32::MIN && ain == bin {
        panic!("gcd({}, {}) overflows", ain, bin);
    }

    let zb = bin.unsigned_abs().trailing_zeros() as usize;
    let mut za = ain.unsigned_abs().trailing_zeros() as usize;
    let mut a = ain.unsigned_abs();
    let mut b = bin.unsigned_abs() >> zb;
    let k = za.min(zb);

    while a != 0 {
        a >>= za;
        let absd = if a >= b { a - b } else { b - a };
        let diff = absd;
        za = diff.trailing_zeros() as usize;
        b = a.min(b);
        a = absd;
    }

    (b << k) as i32
}

/// GCD for i64
pub fn gcd_i64(a: i64, b: i64) -> i64 {
    gcd_bitinteger_i64(a, b)
}

/// Binary GCD for i64 with signed integer handling (matching Julia's algorithm)
pub fn gcd_bitinteger_i64(ain: i64, bin: i64) -> i64 {
    if ain == 0 {
        return bin.abs();
    }
    if bin == 0 {
        return ain.abs();
    }

    if ain == i64::MIN && ain == bin {
        panic!("gcd({}, {}) overflows", ain, bin);
    }

    let zb = bin.unsigned_abs().trailing_zeros() as usize;
    let mut za = ain.unsigned_abs().trailing_zeros() as usize;
    let mut a = ain.unsigned_abs();
    let mut b = bin.unsigned_abs() >> zb;
    let k = za.min(zb);

    while a != 0 {
        a >>= za;
        let absd = if a >= b { a - b } else { b - a };
        let diff = absd;
        za = diff.trailing_zeros() as usize;
        b = a.min(b);
        a = absd;
    }

    (b << k) as i64
}

/// GCD for i128
pub fn gcd_i128(a: i128, b: i128) -> i128 {
    gcd_bitinteger_i128(a, b)
}

/// Binary GCD for i128 with signed integer handling (matching Julia's algorithm)
pub fn gcd_bitinteger_i128(ain: i128, bin: i128) -> i128 {
    if ain == 0 {
        return bin.abs();
    }
    if bin == 0 {
        return ain.abs();
    }

    if ain == i128::MIN && ain == bin {
        panic!("gcd({}, {}) overflows", ain, bin);
    }

    let zb = bin.unsigned_abs().trailing_zeros() as usize;
    let mut za = ain.unsigned_abs().trailing_zeros() as usize;
    let mut a = ain.unsigned_abs();
    let mut b = bin.unsigned_abs() >> zb;
    let k = za.min(zb);

    while a != 0 {
        a >>= za;
        let absd = if a >= b { a - b } else { b - a };
        let diff = absd;
        za = diff.trailing_zeros() as usize;
        b = a.min(b);
        a = absd;
    }

    (b << k) as i128
}

#[cfg(test)]
mod tests {
    use super::*;

    #[test]
    fn test_gcd_u64() {
        assert_eq!(gcd_u64(6, 9), 3);
        assert_eq!(gcd_u64(6, 0), 6);
        assert_eq!(gcd_u64(0, 0), 0);
        assert_eq!(gcd_u64(48, 18), 6);
        assert_eq!(gcd_u64(17, 13), 1);
    }

    #[test]
    fn test_gcd_i64() {
        assert_eq!(gcd_i64(6, 9), 3);
        assert_eq!(gcd_i64(6, -9), 3);
        assert_eq!(gcd_i64(-6, 9), 3);
        assert_eq!(gcd_i64(-6, -9), 3);
        assert_eq!(gcd_i64(6, 0), 6);
        assert_eq!(gcd_i64(0, 0), 0);
    }

    #[test]
    fn test_gcd_u32() {
        assert_eq!(gcd_u32(48, 18), 6);
        assert_eq!(gcd_u32(100, 25), 25);
    }

    #[test]
    fn test_gcd_i32() {
        assert_eq!(gcd_i32(48, -18), 6);
        assert_eq!(gcd_i32(-100, 25), 25);
    }
}

// Function to approximate pi using probability that two numbers are coprime
fn calc_pi(n: u64) -> f64 {
    let mut cnt = 0u64; // Counter for coprime pairs
    // Loop through all pairs (a, b) where 1 <= a, b <= N
    for a in 1..=n {
        for b in 1..=n {
            // Check if a and b are coprime using the fast binary GCD algorithm
            if gcd::gcd_u64(a, b) == 1 {
                cnt += 1;  // Increment counter if coprime
            }
        }
    }
    // Probability that two numbers are coprime
    let prob = cnt as f64 / (n as f64 * n as f64);
    // Approximate pi using the formula: pi ≈ sqrt(6 / prob)
    (6.0 / prob).sqrt()
}

// Main function to run the pi approximation
fn main() {
    let n = 10000u64;   // Number limit for coprimality checking
    let start = std::time::Instant::now();
    let pi = calc_pi(n); // Approximate pi
    let duration = start.elapsed();
    println!("calcPi: {:?}", duration);
    println!("N: {}", n);   // Output N
    println!("pi: {}", pi); // Output approximation of pi
}

rustc -C opt-level=3 gcd.rs && ./gcd
calcPi: 1.2510425s
N: 10000
pi: 3.141534239016629